{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "ename": "ModuleNotFoundError",
     "evalue": "No module named 'svm'",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mModuleNotFoundError\u001b[0m                       Traceback (most recent call last)",
      "Cell \u001b[0;32mIn[1], line 3\u001b[0m\n\u001b[1;32m      1\u001b[0m \u001b[39mimport\u001b[39;00m \u001b[39msys\u001b[39;00m\n\u001b[1;32m      2\u001b[0m \u001b[39mfrom\u001b[39;00m \u001b[39mnumpy\u001b[39;00m \u001b[39mimport\u001b[39;00m \u001b[39m*\u001b[39m\n\u001b[0;32m----> 3\u001b[0m \u001b[39mfrom\u001b[39;00m \u001b[39msvm\u001b[39;00m \u001b[39mimport\u001b[39;00m \u001b[39m*\u001b[39m\n\u001b[1;32m      4\u001b[0m \u001b[39mfrom\u001b[39;00m \u001b[39mos\u001b[39;00m \u001b[39mimport\u001b[39;00m listdir\n\u001b[1;32m      5\u001b[0m \u001b[39mimport\u001b[39;00m \u001b[39mdatetime\u001b[39;00m\n",
      "\u001b[0;31mModuleNotFoundError\u001b[0m: No module named 'svm'"
     ]
    }
   ],
   "source": [
    "import sys\n",
    "from numpy import *\n",
    "from svm import *\n",
    "from os import listdir\n",
    "import datetime"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "class PlattSMO:\n",
    "    def __init__(self, dataMat, classlabels, C, toler, maxIter, **kernelargs):\n",
    "        \"\"\"\n",
    "        构造函数\n",
    "        :param dataMat: 特征向量，二维数组\n",
    "        :param classlabels: 数据标签，一维数组\n",
    "        :param C: 对松弛变量的容忍程度，越大越不容忍\n",
    "        :param toler: 完成一次迭代维度误差要求\n",
    "        :param maxIter: 迭代次数\n",
    "        :param kernelargs: 核参数，特别注意高斯核会有两个参数，但是都是通过这一个变量传进来的\n",
    "        \"\"\"\n",
    "        self.x = array(dataMat)\n",
    "        self.label = array(classlabels).transpose()\n",
    "        self.C = C\n",
    "        self.toler = toler\n",
    "        self.maxIter = maxIter\n",
    "        self.m = shape(dataMat)[0]  # 输入的行数，也即是表示了有多少个输入\n",
    "        self.n = shape(dataMat)[1]  # 列数，表示每一个输入有多少个特征向量\n",
    "        self.alpha = array(zeros(self.m), dtype='float64')  # 初始化alpha\n",
    "        self.b = 0.0\n",
    "        self.eCache = array(zeros((self.m, 2)))  # 错误缓存\n",
    "        self.K = zeros((self.m, self.m), dtype='float64')  # 先求内积，\n",
    "        self.kwargs = kernelargs\n",
    "        self.SV = ()  # 最后保留的支持向量\n",
    "        self.SVIndex = None  # 支持向量的索引\n",
    "        # wx本就是内积本就是输入一个特征向量，和已有的所有特征向量内积（或者说支持向量），然后再加权求和求出一个数\n",
    "        # 训练时候，x即是输入的特征向量，也是所有特征向量的一部分，所以可以这样\n",
    "        # 一个确定xi就是输入x，和所有的特征向量，也就是xj求内积，放在一行当中，那应该是按行取为什么按列\n",
    "        # 为了用矩阵求和，而且其实这个是一个对称矩阵，所以可以直接这样，至于为什么对称，emmm\n",
    "        for i in range(self.m):\n",
    "            for j in range(self.m):\n",
    "                self.K[i, j] = self.kernelTrans(self.x[i, :], self.x[j, :])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def calcEK(self, k):\n",
    "        \"\"\"\n",
    "        计算第k个数据的误差\n",
    "        :param k:\n",
    "        :return:\n",
    "        \"\"\"\n",
    "        # 因为这是训练阶段，用数据集的数据，所以可以直接这样做，这里先把内积全部求了，\n",
    "        fxk = dot(self.alpha * self.label, self.K[:, k]) + self.b\n",
    "        Ek = fxk - float(self.label[k])\n",
    "        return Ek\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    " def updateEK(self, k):\n",
    "        Ek = self.calcEK(k)\n",
    "\n",
    "        self.eCache[k] = [1, Ek]\n",
    "\n",
    "    def selectJ(self, i, Ei):\n",
    "        \"\"\"\n",
    "        在确定了一个参数的前提下，按照最大步长取另一个参数\n",
    "        :param i:\n",
    "        :param Ei:\n",
    "        :return:\n",
    "        \"\"\"\n",
    "        maxE = 0.0\n",
    "        selectJ = 0\n",
    "        Ej = 0.0\n",
    "        validECacheList = nonzero(self.eCache[:, 0])[0]\n",
    "        if len(validECacheList) > 1:\n",
    "            for k in validECacheList:\n",
    "                if k == i: continue\n",
    "                Ek = self.calcEK(k)\n",
    "                deltaE = abs(Ei - Ek)\n",
    "                if deltaE > maxE:\n",
    "                    selectJ = k\n",
    "                    maxE = deltaE\n",
    "                    Ej = Ek\n",
    "            return selectJ, Ej\n",
    "        else:\n",
    "            selectJ = selectJrand(i, self.m)\n",
    "            Ej = self.calcEK(selectJ)\n",
    "            return selectJ, Ej"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def innerL(self, i):\n",
    "        \"\"\"\n",
    "        选择参数之后，更新参数\n",
    "        \"\"\"\n",
    "        Ei = self.calcEK(i)\n",
    "\n",
    "        if (self.label[i] * Ei < -self.toler and self.alpha[i] < self.C) or \\\n",
    "                (self.label[i] * Ei > self.toler and self.alpha[i] > 0):\n",
    "\n",
    "            self.updateEK(i)\n",
    "\n",
    "            j, Ej = self.selectJ(i, Ei)\n",
    "\n",
    "            alphaIOld = self.alpha[i].copy()\n",
    "            alphaJOld = self.alpha[j].copy()\n",
    "\n",
    "            if self.label[i] != self.label[j]:\n",
    "                L = max(0, self.alpha[j] - self.alpha[i])\n",
    "                H = min(self.C, self.C + self.alpha[j] - self.alpha[i])\n",
    "            else:\n",
    "                L = max(0, self.alpha[j] + self.alpha[i] - self.C)\n",
    "                H = min(self.C, self.alpha[i] + self.alpha[j])\n",
    "            if L == H:\n",
    "                return 0\n",
    "\n",
    "            eta = 2 * self.K[i, j] - self.K[i, i] - self.K[j, j]\n",
    "            if eta >= 0:\n",
    "                return 0\n",
    "\n",
    "            self.alpha[j] -= self.label[j] * (Ei - Ej) / eta\n",
    "            self.alpha[j] = clipAlpha(self.alpha[j], H, L)\n",
    "            self.updateEK(j)\n",
    "\n",
    "            if abs(alphaJOld - self.alpha[j]) < 0.00001:  # 目标迭代完成，不需要再迭代，而且所有参数已经保留，理论上是不应该保留当前的参数的，但是也正因为反正相差不多，所以可以\n",
    "                return 0\n",
    "\n",
    "            self.alpha[i] += self.label[i] * self.label[j] * (alphaJOld - self.alpha[j])\n",
    "            self.updateEK(i)\n",
    "\n",
    "            b1 = self.b - Ei - self.label[i] * self.K[i, i] * (self.alpha[i] - alphaIOld) - \\\n",
    "                 self.label[j] * self.K[i, j] * (self.alpha[j] - alphaJOld)\n",
    "            b2 = self.b - Ej - self.label[i] * self.K[i, j] * (self.alpha[i] - alphaIOld) - \\\n",
    "                 self.label[j] * self.K[j, j] * (self.alpha[j] - alphaJOld)\n",
    "            if 0 < self.alpha[i] and self.alpha[i] < self.C:\n",
    "                self.b = b1\n",
    "            elif 0 < self.alpha[j] and self.alpha[j] < self.C:\n",
    "                self.b = b2\n",
    "            else:\n",
    "                self.b = (b1 + b2) / 2.0\n",
    "            return 1\n",
    "        else:\n",
    "            return 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    " def smoP(self):\n",
    "        \"\"\"\n",
    "        外层大循环，\n",
    "        :return:\n",
    "        \"\"\"\n",
    "        iter = 0\n",
    "        entrySet = True\n",
    "        alphaPairChanged = 0\n",
    "        while iter < self.maxIter and ((alphaPairChanged > 0) or (entrySet)):\n",
    "            alphaPairChanged = 0\n",
    "            if entrySet:\n",
    "                for i in range(self.m):\n",
    "                    alphaPairChanged += self.innerL(i)\n",
    "                iter += 1\n",
    "            else:\n",
    "                nonBounds = nonzero((self.alpha > 0) * (self.alpha < self.C))[0]\n",
    "                for i in nonBounds:\n",
    "                    alphaPairChanged += self.innerL(i)\n",
    "                iter += 1\n",
    "            if entrySet:\n",
    "                entrySet = False\n",
    "            elif alphaPairChanged == 0:\n",
    "                entrySet = True\n",
    "        # 保存模型参数\n",
    "        self.SVIndex = nonzero(self.alpha)[0]  # 取非0部分把应该\n",
    "        self.SV = self.x[self.SVIndex]\n",
    "        self.SVAlpha = self.alpha[self.SVIndex]\n",
    "        self.SVLabel = self.label[self.SVIndex]\n",
    "\n",
    "        # 清空中间变量\n",
    "        self.x = None\n",
    "        self.K = None\n",
    "        self.label = None\n",
    "        self.alpha = None\n",
    "        self.eCache = None\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#   def K(self,i,j):\n",
    "    #       return self.x[i,:]*self.x[j,:].T\n",
    "    def kernelTrans(self, x, z):\n",
    "        \"\"\"\n",
    "        核函数说到底就是求内积，求输入x和已有数据标签的内积，最后的结果是一个常数\n",
    "        内积有两个，一个是向量形式的对应相乘相加一个是矩阵形式的矩阵乘积\n",
    "        输入是两个要求内积的一维数组，也就是说要拆成一个个来做\n",
    "        :param x:\n",
    "        :param z:\n",
    "        :return:\n",
    "        \"\"\"\n",
    "        if array(x).ndim != 1 or array(x).ndim != 1:\n",
    "            raise Exception(\"input vector is not 1 dim\")\n",
    "        if self.kwargs['name'] == 'linear':\n",
    "            return sum(x * z)\n",
    "        elif self.kwargs['name'] == 'rbf':\n",
    "            theta = self.kwargs['theta']\n",
    "            return exp(sum((x - z) * (x - z)) / (-1 * theta ** 2))\n",
    "\n",
    "    def calcw(self):\n",
    "        \"\"\"\n",
    "        计算w，结果是一个数组，长度是特征向量的长度\n",
    "        :return:\n",
    "        \"\"\"\n",
    "        for i in range(self.m):\n",
    "            self.w += dot(self.alpha[i] * self.label[i], self.x[i, :])\n",
    "\n",
    "    def predict(self, testData):\n",
    "        \"\"\"\n",
    "        输入待预测的数据，输出结果\n",
    "        :param testData: 待预测数据，要是数组形式，二维数组\n",
    "        :return: 一个列表，包含结果\n",
    "        \"\"\"\n",
    "        test = array(testData)\n",
    "        # return (test * self.w + self.b).getA()\n",
    "        result = []\n",
    "        m = shape(test)[0]\n",
    "        for i in range(m):\n",
    "            tmp = self.b\n",
    "            for j in range(len(self.SVIndex)):\n",
    "                # wx+b,w和SVAlpha，Label以及核函数有关\n",
    "                # 求和可以是直接求和，也可以转成矩阵，这里就是直接求和\n",
    "                # 计算支持向量的数目，用他们来进行估计\n",
    "                tmp += self.SVAlpha[j] * self.SVLabel[j] * self.kernelTrans(self.SV[j], test[i, :])\n",
    "\n",
    "            # 不可分的情况下，其实就可以，直接指定想要的情况，如果需要的话，工程中会有这个要求\n",
    "            while tmp == 0:\n",
    "                tmp = random.uniform(-1, 1)\n",
    "            if tmp > 0:\n",
    "                tmp = 1\n",
    "            else:\n",
    "                tmp = -1\n",
    "            result.append(tmp)\n",
    "        return result"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plotBestfit(data, label, w, b):\n",
    "    import matplotlib.pyplot as plt\n",
    "    n = shape(data)[0]\n",
    "    fig = plt.figure()\n",
    "    ax = fig.add_subplot(111)\n",
    "    x1 = []\n",
    "    x2 = []\n",
    "    y1 = []\n",
    "    y2 = []\n",
    "    for i in range(n):\n",
    "        if int(label[i]) == 1:\n",
    "            x1.append(data[i][0])\n",
    "            y1.append(data[i][1])\n",
    "        else:\n",
    "            x2.append(data[i][0])\n",
    "            y2.append(data[i][1])\n",
    "    ax.scatter(x1, y1, s=10, c='red', marker='s')\n",
    "    ax.scatter(x2, y2, s=10, c='green', marker='s')\n",
    "    x = arange(-2, 10, 0.1)\n",
    "    y = ((-b - w[0] * x) / w[1])\n",
    "    plt.plot(x, y)\n",
    "    plt.xlabel('X')\n",
    "    plt.ylabel('y')\n",
    "    plt.show()\n",
    "\n",
    "\n",
    "def loadImage(dir, maps=None):\n",
    "    dirList = listdir(dir)\n",
    "    data = []\n",
    "    label = []\n",
    "    for file in dirList:\n",
    "        label.append(file.split('_')[0])\n",
    "        lines = open(dir + '/' + file).readlines()\n",
    "        row = len(lines)\n",
    "        col = len(lines[0].strip())\n",
    "        line = []\n",
    "        for i in range(row):\n",
    "            for j in range(col):\n",
    "                line.append(float(lines[i][j]))\n",
    "        data.append(line)\n",
    "        if maps != None:\n",
    "            label[-1] = float(maps[label[-1]])\n",
    "        else:\n",
    "            label[-1] = float(label[-1])\n",
    "    return array(data), array(label)\n",
    "\n",
    "\n",
    "def main():\n",
    "    '''\n",
    "    data,label = loadDataSet('testSetRBF.txt')\n",
    "    smo = PlattSMO(data,label,200,0.0001,10000,name = 'rbf',theta = 1.3)\n",
    "    smo.smoP()\n",
    "    smo.calcw()\n",
    "    print smo.predict(data)\n",
    "    '''\n",
    "    # 数据读取\n",
    "    # 特征是二维数组，标签是一维数组，记住这个就可以，而且标签不是实际标签，而都要转变成-1,1，正类和负类。，谁正谁负没关系。\n",
    "    maps = {'1': 1.0, '9': -1.0}\n",
    "    data, label = loadImage(\"digits/trainingDigits\", maps)\n",
    "    test, testLabel = loadImage(\"digits/testDigits\", maps)\n",
    "\n",
    "    # 训练\n",
    "    smo = PlattSMO(data, label, 1, 0.0001, 10000, name='rbf', theta=20)\n",
    "    begin = datetime.datetime.now()\n",
    "    smo.smoP()\n",
    "    end = datetime.datetime.now()\n",
    "    time_sub = end - begin\n",
    "    print(\"keneral_svm time\", time_sub.total_seconds())\n",
    "    print(len(smo.SVIndex))\n",
    "\n",
    "    # 预测\n",
    "    testResult = smo.predict(test)\n",
    "    m = shape(test)[0]\n",
    "    count = 0.0\n",
    "    for i in range(m):\n",
    "        if testLabel[i] != testResult[i]:\n",
    "            count += 1\n",
    "    print(\"classfied error rate is:\", count / m)\n",
    "    # smo.kernelTrans(data,smo.SV[0])\n",
    "\n",
    "\n",
    "if __name__ == \"__main__\":\n",
    "    sys.exit(main())"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.6"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
